library(Seurat)
library(Signac)
library(tidyverse)
library(ggplot2)
library(ggpubr)
library(magick)
library(knitr)
library(kableExtra)
library(devtools)
library(harmony)
library(patchwork)
library(kableExtra)
# Paths
path_to_obj <- ("~/Documents/multiome_tonsil_Lucia/results/R_objects/11.tonsil_multiome_integrated_without_doublets_normalized.rds")
path_to_save <- ("~/Documents/multiome_tonsil_Lucia/results/R_objects/12.tonsil_multiome_bcells_without_doublets_normalized.rds")
tonsil_wnn_bcell <- readRDS(path_to_obj)
vars <- str_subset(colnames(tonsil_wnn_bcell@meta.data), "^wsnn_res")
clusters_gg <- purrr::map(vars, function(x) {
p <- DimPlot(
tonsil_wnn_bcell,
group.by = x,
reduction = "wnn.umap",
pt.size = 0.1, label = T
)
p
})
clusters_gg
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
tonsil_wnn_bcell$is_tcell <-
tonsil_wnn_bcell$wsnn_res.0.01 == "1"
tonsil_wnn_bcell <- subset(tonsil_wnn_bcell, subset = is_tcell == FALSE)
DimPlot(tonsil_wnn_bcell, reduction = "wnn.umap", group.by = "wsnn_res.0.01", label = TRUE, pt.size = 0.1)
We exclude the first dimension as this is typically correlated with sequencing depth Cells cluster completely separately in ATAC without harmony; so run harmony after SVD
RunSVD LSI
DefaultAssay(tonsil_wnn_bcell) <- "peaks"
tonsil_wnn_bcell <- RunTFIDF(tonsil_wnn_bcell)
## Performing TF-IDF normalization
tonsil_wnn_bcell <- FindTopFeatures(tonsil_wnn_bcell, min.cutoff = "q0")
tonsil_wnn_bcell <- RunSVD(tonsil_wnn_bcell)
## Running SVD
## Scaling cell embeddings
Compute the correlation between total counts and each reduced dimension component.
LSI component is typically highly correlated with sequencing depth. The first LSI component often captures sequencing depth (technical variation) rather than biological variation. If this is the case, the component should be removed from downstream analysis. We can assess the correlation between each LSI component and sequencing depth using the DepthCor() function:
For scRNA-seq data we don’t typically observe such a strong relationship between the first PC and sequencing depth, and so usually retain the first PC in downstream analyses.
DepthCor(tonsil_wnn_bcell)
Here we see there is a very strong correlation between the first LSI component and the total number of counts for the cell, so we will perform downstream steps without this component.
tonsil_wnn_bcell <- RunUMAP(
tonsil_wnn_bcell,
dims = 2:40,
reduction = "lsi",
reduction.name = "umap.atac",
reduction.key = "atacUMAP_"
)
## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per session
## 11:55:55 UMAP embedding parameters a = 0.9922 b = 1.112
## 11:55:56 Read 46120 rows and found 39 numeric columns
## 11:55:56 Using Annoy for neighbor search, n_neighbors = 30
## 11:55:57 Building Annoy index with metric = cosine, n_trees = 50
## 0% 10 20 30 40 50 60 70 80 90 100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 11:56:14 Writing NN index file to temp file /var/folders/kz/y10np0cj213fybqz181z9lghbzr42p/T//RtmpP3rndT/file3fc26f8884bd
## 11:56:15 Searching Annoy index using 1 thread, search_k = 3000
## 11:56:54 Annoy recall = 100%
## 11:57:22 Commencing smooth kNN distance calibration using 1 thread
## 11:57:27 Initializing from normalized Laplacian + noise
## 11:57:31 Commencing optimization for 200 epochs, with 1988358 positive edges
## 11:58:04 Optimization finished
atac.umap<-DimPlot(
tonsil_wnn_bcell,
reduction = "umap.atac",
group.by = "library_name",
pt.size = 0.1
) + ggtitle('scATAC UMAP') + NoLegend()
atac.umap
#split_by: library ,edad, genero
DefaultAssay(tonsil_wnn_bcell) <- "RNA"
tonsil_wnn_bcell <- NormalizeData(
tonsil_wnn_bcell,
normalization.method = "LogNormalize",
scale.factor = 1e4
)
tonsil_wnn_bcell <- tonsil_wnn_bcell %>%
FindVariableFeatures(nfeatures = 3000) %>%
ScaleData() %>%
RunPCA()
PCAPlot(tonsil_wnn_bcell,
group.by = "library_name")
ElbowPlot(object = tonsil_wnn_bcell)
find variable genes
tonsil_wnn_bcell <- RunUMAP(
tonsil_wnn_bcell,
dims = 1:30,
reduction = "pca",
reduction.name = "umap.rna",
reduction.key = "rnaUMAP_"
)
## 12:00:32 UMAP embedding parameters a = 0.9922 b = 1.112
## 12:00:32 Read 46120 rows and found 30 numeric columns
## 12:00:32 Using Annoy for neighbor search, n_neighbors = 30
## 12:00:32 Building Annoy index with metric = cosine, n_trees = 50
## 0% 10 20 30 40 50 60 70 80 90 100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 12:00:36 Writing NN index file to temp file /var/folders/kz/y10np0cj213fybqz181z9lghbzr42p/T//RtmpP3rndT/file3fc2448a7bba
## 12:00:36 Searching Annoy index using 1 thread, search_k = 3000
## 12:00:51 Annoy recall = 100%
## 12:00:52 Commencing smooth kNN distance calibration using 1 thread
## 12:00:55 Initializing from normalized Laplacian + noise
## 12:00:58 Commencing optimization for 200 epochs, with 2043478 positive edges
## 12:01:28 Optimization finished
rna.umap<-DimPlot(
tonsil_wnn_bcell,
reduction = "umap.rna",
group.by = "library_name",
pt.size = 0.1) + NoLegend() + ggtitle('scRNA UMAP')
rna.umap
atac.umap + rna.umap
hacer primero harmony, quitar batch effect. atac y rna. harmony
Pass the Seurat object to the RunHarmony function and specify which variable to integrate out. A Seurat object is returned with corrected Harmony coordinates.
DefaultAssay(tonsil_wnn_bcell) <- "peaks"
tonsil_wnn_bcell <- RunHarmony(
object = tonsil_wnn_bcell,
reduction = "lsi",
dims = 2:40,
group.by.vars = "library_name",
assay.use = "peaks",
project.dim = FALSE,
reduction.save = "harmony_peaks"
)
## Harmony 1/10
## Harmony 2/10
## Harmony 3/10
## Harmony 4/10
## Harmony 5/10
## Harmony 6/10
## Harmony 7/10
## Harmony 8/10
## Harmony converged after 8 iterations
tonsil_wnn_bcell <- RunUMAP(
tonsil_wnn_bcell,
dims = 2:40,
reduction = "harmony_peaks",
reduction.name = "umap.atac",
reduction.key = "atacUMAP_"
)
## 12:04:55 UMAP embedding parameters a = 0.9922 b = 1.112
## 12:04:55 Read 46120 rows and found 39 numeric columns
## 12:04:55 Using Annoy for neighbor search, n_neighbors = 30
## 12:04:55 Building Annoy index with metric = cosine, n_trees = 50
## 0% 10 20 30 40 50 60 70 80 90 100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 12:05:00 Writing NN index file to temp file /var/folders/kz/y10np0cj213fybqz181z9lghbzr42p/T//RtmpP3rndT/file3fc27a76b6ea
## 12:05:00 Searching Annoy index using 1 thread, search_k = 3000
## 12:05:17 Annoy recall = 100%
## 12:05:18 Commencing smooth kNN distance calibration using 1 thread
## 12:05:20 Initializing from normalized Laplacian + noise
## 12:05:22 Commencing optimization for 200 epochs, with 2053750 positive edges
## 12:05:49 Optimization finished
Harm_peak<-DimPlot(
tonsil_wnn_bcell,
reduction = "umap.atac",
group.by = "library_name",
pt.size = 0.1
) + NoLegend() + ggtitle('Peak Harmony')
tonsil_wnn_bcell <- RunUMAP(
tonsil_wnn_bcell,
dims = 2:40,
reduction = "harmony_rna",
reduction.name = "umap.rna",
reduction.key = "rnaUMAP_"
)
## 12:07:14 UMAP embedding parameters a = 0.9922 b = 1.112
## 12:07:14 Read 46120 rows and found 39 numeric columns
## 12:07:14 Using Annoy for neighbor search, n_neighbors = 30
## 12:07:14 Building Annoy index with metric = cosine, n_trees = 50
## 0% 10 20 30 40 50 60 70 80 90 100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 12:07:18 Writing NN index file to temp file /var/folders/kz/y10np0cj213fybqz181z9lghbzr42p/T//RtmpP3rndT/file3fc26f5508a2
## 12:07:18 Searching Annoy index using 1 thread, search_k = 3000
## 12:07:34 Annoy recall = 100%
## 12:07:34 Commencing smooth kNN distance calibration using 1 thread
## 12:07:37 Initializing from normalized Laplacian + noise
## 12:07:39 Commencing optimization for 200 epochs, with 2100454 positive edges
## 12:08:07 Optimization finished
Harm_rna<-DimPlot(
tonsil_wnn_bcell,
reduction = "umap.rna",
group.by = "library_name",
pt.size = 0.1
) + NoLegend() + ggtitle('RNA Harmony')
Harm_peak+Harm_rna
FindNeighbors
Constructs a Shared Nearest Neighbor (SNN) Graph for a given dataset. We first determine the k-nearest neighbors of each cell. We use this knn graph to construct the SNN graph by calculating the neighborhood overlap (Jaccard index) between every cell and its k.param nearest neighbors.
# build a joint neighbor graph using both assays
tonsil_wnn_bcell <- FindMultiModalNeighbors(
object = tonsil_wnn_bcell,
reduction.list = list("harmony_peaks", "harmony_rna"),
dims.list = list(2:40, 1:30), modality.weight.name = "Joint_snn_umap"
)
## Calculating cell-specific modality weights
## Finding 20 nearest neighbors for each modality.
## Calculating kernel bandwidths
## Warning in FindMultiModalNeighbors(object = tonsil_wnn_bcell, reduction.list
## = list("harmony_peaks", : The number of provided modality.weight.name is not
## equal to the number of modalities. peaks.weight RNA.weight are used to store the
## modality weights
## Finding multimodal neighbors
## Constructing multimodal KNN graph
## Constructing multimodal SNN graph
# build a joint UMAP visualization
tonsil_wnn_bcell <- RunUMAP(
object = tonsil_wnn_bcell,
nn.name = "weighted.nn",
reduction.name = "wnn.umap",
reduction.key = "wnnUMAP_")
## 12:10:39 UMAP embedding parameters a = 0.9922 b = 1.112
## 12:10:41 Commencing smooth kNN distance calibration using 1 thread
## 12:10:43 Initializing from normalized Laplacian + noise
## 12:10:45 Commencing optimization for 200 epochs, with 1475168 positive edges
## 12:11:15 Optimization finished
joint.umap<- DimPlot(tonsil_wnn_bcell, label = FALSE, group.by = "library_name", pt.size = 0.1, reduction = "wnn.umap") + plot_annotation(title = 'Joint UMAP')+ ggtitle('Joint UMAP by library name') + NoLegend()
joint.umap
joint.umap_age<- DimPlot(tonsil_wnn_bcell, label = FALSE, split.by = "age_group", pt.size = 0.1, reduction = "wnn.umap") + plot_annotation(title = 'Joint UMAP ')+ ggtitle('Joint UMAP age ')
joint.umap_hospital<- DimPlot(tonsil_wnn_bcell, label = FALSE, split.by = "hospital", pt.size = 0.1, reduction = "wnn.umap") + plot_annotation(title = 'Joint UMAP ')+ ggtitle('Joint UMAP age ')
joint.umap_age
joint.umap_hospital
#find cluster algorithm 3 = SLM algorithm
tonsil_wnn_bcell <- FindClusters(tonsil_wnn_bcell, resolution = c(0.005,0.01,0.05,0.075,0.1,0.25,0.5,0.75),algorithm = 3, graph.name = "wsnn",verbose = FALSE)
print(colnames(tonsil_wnn_bcell@meta.data))
## [1] "lib_name_barcode" "orig.ident" "nCount_RNA"
## [4] "nFeature_RNA" "nCount_ATAC" "nFeature_ATAC"
## [7] "nucleosome_signal" "nucleosome_percentile" "TSS.enrichment"
## [10] "TSS.percentile" "tss.level" "percent.mt"
## [13] "percent_ribo" "nCount_peaks" "nFeature_peaks"
## [16] "library_name" "donor_id" "sex"
## [19] "age" "age_group" "hospital"
## [22] "assay" "barcodes" "doublet_scores"
## [25] "predicted_doublets" "peaks.weight" "RNA.weight"
## [28] "wsnn_res.0.005" "wsnn_res.0.01" "seurat_clusters"
## [31] "sub.cluster_0.25" "sub.cluster0_0.5" "is_doublet"
## [34] "wsnn_res.0.05" "wsnn_res.0.75" "wsnn_res.0.075"
## [37] "is_tcell" "wsnn_res.0.1" "wsnn_res.0.25"
## [40] "wsnn_res.0.5"
vars <- str_subset(colnames(tonsil_wnn_bcell@meta.data), "^wsnn_res")
clusters_gg <- purrr::map(vars, function(x) {
p <- DimPlot(
tonsil_wnn_bcell,
group.by = x,
reduction = "wnn.umap",
pt.size = 0.1, label = T
)
p
})
clusters_gg
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
##
## [[6]]
##
## [[7]]
##
## [[8]]
Idents(tonsil_wnn_bcell)<-"wsnn_res.0.05"
tonsil_markers_05<-FindAllMarkers(object = tonsil_wnn_bcell, only.pos = TRUE)
## Calculating cluster 0
## Calculating cluster 1
## Calculating cluster 2
## Calculating cluster 3
## Calculating cluster 4
## Calculating cluster 5
## Calculating cluster 6
write.csv(tonsil_markers_05,file=paste0("~/Documents/multiome_tonsil_Lucia/results/tables/", "tonsil_markers_bcell_05.csv"))
Resolution 0.01
tonsil_markers_05 %>% group_by(cluster) %>% top_n(n = 10, wt = avg_log2FC) %>% write.csv(.,file=paste0("~/Documents/multiome_tonsil_Lucia/results/tables/", "top10_tonsil_markers_bcell_05.csv"))
tonsil_markers_05 %>% group_by(cluster) %>% top_n(n = 5, wt = avg_log2FC) %>% write.csv(.,file=paste0("~/Documents/multiome_tonsil_Lucia/results/tables/", "top5_tonsil_markers_bcell_05.csv"))
top5_tonsil_markers_05<-tonsil_markers_05 %>% group_by(cluster) %>% top_n(n = 5, wt = avg_log2FC)
top10_tonsil_markers_05<-tonsil_markers_05 %>% group_by(cluster) %>% top_n(n = 10, wt = avg_log2FC)
df_top5<-as.data.frame(top5_tonsil_markers_05)
kbl(df_top5,caption = "Table of the top 5 marker of each cluster resolution 0.005") %>%
kable_paper("striped", full_width = F)
| p_val | avg_log2FC | pct.1 | pct.2 | p_val_adj | cluster | gene |
|---|---|---|---|---|---|---|
| 0 | 2.439037 | 0.994 | 0.768 | 0 | 0 | BANK1 |
| 0 | 2.310305 | 0.618 | 0.170 | 0 | 0 | COL19A1 |
| 0 | 2.002782 | 0.806 | 0.357 | 0 | 0 | CAMK1D |
| 0 | 1.857727 | 0.729 | 0.316 | 0 | 0 | MAML2 |
| 0 | 1.791176 | 0.639 | 0.227 | 0 | 0 | PARP15 |
| 0 | 2.918946 | 0.955 | 0.166 | 0 | 1 | HMGB2 |
| 0 | 2.798930 | 0.967 | 0.254 | 0 | 1 | TUBA1B |
| 0 | 2.694975 | 0.967 | 0.262 | 0 | 1 | H2AFZ |
| 0 | 2.422139 | 0.851 | 0.335 | 0 | 1 | HIST1H4C |
| 0 | 2.412212 | 0.815 | 0.043 | 0 | 1 | TOP2A |
| 0 | 2.612691 | 0.835 | 0.260 | 0 | 2 | MAML3 |
| 0 | 2.542482 | 0.983 | 0.299 | 0 | 2 | AC023590.1 |
| 0 | 2.459488 | 0.819 | 0.165 | 0 | 2 | AC104170.1 |
| 0 | 2.257369 | 0.929 | 0.250 | 0 | 2 | RAPGEF5 |
| 0 | 2.156023 | 0.796 | 0.249 | 0 | 2 | LHFPL2 |
| 0 | 2.505368 | 0.907 | 0.042 | 0 | 3 | FYB1 |
| 0 | 2.485483 | 0.816 | 0.038 | 0 | 3 | INPP4B |
| 0 | 2.066960 | 0.726 | 0.010 | 0 | 3 | THEMIS |
| 0 | 1.998146 | 0.902 | 0.147 | 0 | 3 | PRKCH |
| 0 | 1.978889 | 0.807 | 0.010 | 0 | 3 | BCL11B |
| 0 | 5.751547 | 0.498 | 0.151 | 0 | 4 | IGHGP |
| 0 | 5.587754 | 0.697 | 0.326 | 0 | 4 | IGHG3 |
| 0 | 5.666937 | 0.960 | 0.918 | 0 | 4 | IGKC |
| 0 | 5.984929 | 0.648 | 0.439 | 0 | 4 | IGHA1 |
| 0 | 5.614213 | 0.896 | 0.762 | 0 | 4 | IGLC2 |
| 0 | 4.211622 | 0.735 | 0.056 | 0 | 5 | SLC8A1 |
| 0 | 3.942508 | 0.687 | 0.009 | 0 | 5 | LYZ |
| 0 | 3.486982 | 0.687 | 0.086 | 0 | 5 | AOAH |
| 0 | 3.469585 | 0.197 | 0.014 | 0 | 5 | SPRR3 |
| 0 | 3.357581 | 0.695 | 0.004 | 0 | 5 | PLXDC2 |
| 0 | 4.279665 | 0.976 | 0.081 | 0 | 6 | LINC01374 |
| 0 | 3.938582 | 0.920 | 0.011 | 0 | 6 | LINC01478 |
| 0 | 3.801066 | 0.992 | 0.142 | 0 | 6 | RUNX2 |
| 0 | 3.708258 | 0.932 | 0.009 | 0 | 6 | FAM160A1 |
| 0 | 3.361241 | 0.880 | 0.004 | 0 | 6 | DACH1 |
markers_gg <- function(x){purrr::map(x, function(x) {
p <- FeaturePlot(
tonsil_wnn_bcell,
features = x,
reduction = "wnn.umap",
pt.size = 0.1
)
p
})}
MARKERS
Immature B cells express CD19, CD 20, CD34, CD38, and CD45R, T-cell receptor/CD3 complex (TCR/CD3 complex)
DZ: SUGCT, CXCR4, AICDA
LZ: CD83, BCL2A1
GC total: MEF2B, BCL6, IRF4
PC: PRDM1, SLAMF7, MZB1, FKBP11
monocytes_markers<-c("LYZ","S100A8")
naive_markers<-c("CD79A", "CD79B", "BLNK")
bib_Bcell_markers<-c("CD19","CR2","MS4A1","RALGPS2","CD79A")
bib_Tcell_markers<-c("CD3E","CD4","CD8A","FOXP3","IL17A")
markers_bcell<-c("BANK1","ARHGAP24","ADAM28","MARCH1","RAPGEF5","AFF2","RGS13","LPP","IGHG1","IGLC1","SLC8A1","LYZ","PLXDC2","FAM160A1","IGHA1","IGLC2", "SETBP1","ENTPD1","COL19A1","CCSER1")
markers_tcell<-c("INPP4B","FYB1","LEF1","IL7R","IL6ST","CCL5","GNLY","NKG7","DTHD1","RUNX2", "FOXP3","CD8A","IL17A","CD2")
CD8+ T cell markers:“CD3D”, “CD8A” NK cell markers:“GNLY”, “NKG7”
markers_gg (bib_Tcell_markers)
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
markers_gg(bib_Bcell_markers)
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
markers_gg(monocytes_markers)
## [[1]]
##
## [[2]]
markers_gg (naive_markers)
## [[1]]
##
## [[2]]
##
## [[3]]
m<-c("PRDM1","XBP1","IRF4","MEF2B","BCL6")
DZ<-c("SUGCT", "CXCR4", "AICDA")
LZ<- c("CD83","BCL2A1")
GC<- c("MEF2B", "BCL6","IRF4")
PC<- c("PRDM1","SLAMF7", "MZB1", "FKBP11")
markers_gg(m)
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
markers_gg(DZ)
## [[1]]
##
## [[2]]
##
## [[3]]
markers_gg(LZ)
## [[1]]
##
## [[2]]
markers_gg(GC)
## [[1]]
##
## [[2]]
##
## [[3]]
markers_gg(PC)
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
saveRDS(tonsil_wnn_bcell, path_to_save)
sessionInfo()
## R version 4.1.2 (2021-11-01)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
##
## locale:
## [1] es_ES.UTF-8/es_ES.UTF-8/es_ES.UTF-8/C/es_ES.UTF-8/es_ES.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] patchwork_1.1.1 harmony_0.1.0 Rcpp_1.0.8 devtools_2.4.3
## [5] usethis_2.1.5 kableExtra_1.3.4 knitr_1.36 magick_2.7.3
## [9] ggpubr_0.4.0 forcats_0.5.1 stringr_1.4.0 dplyr_1.0.7
## [13] purrr_0.3.4 readr_2.1.1 tidyr_1.1.4 tibble_3.1.6
## [17] ggplot2_3.3.5 tidyverse_1.3.1 Signac_1.5.0 SeuratObject_4.0.4
## [21] Seurat_4.0.6 BiocStyle_2.22.0
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.2 reticulate_1.22 tidyselect_1.1.1
## [4] htmlwidgets_1.5.4 grid_4.1.2 docopt_0.7.1
## [7] BiocParallel_1.28.3 Rtsne_0.15 munsell_0.5.0
## [10] codetools_0.2-18 ica_1.0-2 future_1.23.0
## [13] miniUI_0.1.1.1 withr_2.4.3 colorspace_2.0-2
## [16] highr_0.9 rstudioapi_0.13 stats4_4.1.2
## [19] ROCR_1.0-11 ggsignif_0.6.3 tensor_1.5
## [22] listenv_0.8.0 labeling_0.4.2 slam_0.1-49
## [25] GenomeInfoDbData_1.2.7 polyclip_1.10-0 farver_2.1.0
## [28] rprojroot_2.0.2 parallelly_1.30.0 vctrs_0.3.8
## [31] generics_0.1.1 xfun_0.29 lsa_0.73.2
## [34] ggseqlogo_0.1 R6_2.5.1 GenomeInfoDb_1.30.0
## [37] cachem_1.0.6 bitops_1.0-7 spatstat.utils_2.3-0
## [40] assertthat_0.2.1 promises_1.2.0.1 scales_1.1.1
## [43] gtable_0.3.0 globals_0.14.0 processx_3.5.2
## [46] goftest_1.2-3 rlang_0.4.12 systemfonts_1.0.3
## [49] RcppRoll_0.3.0 splines_4.1.2 rstatix_0.7.0
## [52] lazyeval_0.2.2 spatstat.geom_2.3-1 broom_0.7.10
## [55] BiocManager_1.30.16 yaml_2.2.1 reshape2_1.4.4
## [58] abind_1.4-5 modelr_0.1.8 backports_1.4.1
## [61] httpuv_1.6.4 tools_4.1.2 bookdown_0.24
## [64] ellipsis_0.3.2 spatstat.core_2.3-2 jquerylib_0.1.4
## [67] RColorBrewer_1.1-2 BiocGenerics_0.40.0 sessioninfo_1.2.2
## [70] ggridges_0.5.3 plyr_1.8.6 zlibbioc_1.40.0
## [73] RCurl_1.98-1.5 prettyunits_1.1.1 ps_1.6.0
## [76] rpart_4.1-15 deldir_1.0-6 pbapply_1.5-0
## [79] cowplot_1.1.1 S4Vectors_0.32.3 zoo_1.8-9
## [82] haven_2.4.3 ggrepel_0.9.1 cluster_2.1.2
## [85] fs_1.5.2 magrittr_2.0.1 RSpectra_0.16-0
## [88] data.table_1.14.2 scattermore_0.7 lmtest_0.9-39
## [91] reprex_2.0.1 RANN_2.6.1 SnowballC_0.7.0
## [94] fitdistrplus_1.1-6 matrixStats_0.61.0 pkgload_1.2.4
## [97] hms_1.1.1 mime_0.12 evaluate_0.14
## [100] xtable_1.8-4 sparsesvd_0.2 readxl_1.3.1
## [103] IRanges_2.28.0 gridExtra_2.3 testthat_3.1.1
## [106] compiler_4.1.2 KernSmooth_2.23-20 crayon_1.4.2
## [109] htmltools_0.5.2 mgcv_1.8-38 later_1.3.0
## [112] tzdb_0.2.0 lubridate_1.8.0 DBI_1.1.2
## [115] tweenr_1.0.2 dbplyr_2.1.1 MASS_7.3-54
## [118] Matrix_1.3-4 car_3.0-12 cli_3.1.1
## [121] parallel_4.1.2 igraph_1.2.10 GenomicRanges_1.46.1
## [124] pkgconfig_2.0.3 plotly_4.10.0 spatstat.sparse_2.1-0
## [127] xml2_1.3.3 svglite_2.0.0 bslib_0.3.1
## [130] webshot_0.5.2 XVector_0.34.0 rvest_1.0.2
## [133] callr_3.7.0 digest_0.6.29 sctransform_0.3.2
## [136] RcppAnnoy_0.0.19 spatstat.data_2.1-2 Biostrings_2.62.0
## [139] rmarkdown_2.11 cellranger_1.1.0 leiden_0.3.9
## [142] fastmatch_1.1-3 uwot_0.1.11 shiny_1.7.1
## [145] Rsamtools_2.10.0 lifecycle_1.0.1 nlme_3.1-153
## [148] jsonlite_1.7.3 carData_3.0-4 limma_3.50.0
## [151] desc_1.4.0 viridisLite_0.4.0 fansi_1.0.2
## [154] pillar_1.6.5 lattice_0.20-45 pkgbuild_1.3.0
## [157] fastmap_1.1.0 httr_1.4.2 survival_3.2-13
## [160] remotes_2.4.2 glue_1.6.1 qlcMatrix_0.9.7
## [163] png_0.1-7 ggforce_0.3.3 stringi_1.7.6
## [166] sass_0.4.0 memoise_2.0.1 irlba_2.3.5
## [169] future.apply_1.8.1